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We discuss two problems which are particularly challenging for approximations in time-dependent 
density functional theory (TDDFT) to capture: momentum-distributions in ionization processes, 
and memory-dependence in real-time dynamics. We propose an extension of TDDFT to phase- 
space densities, discuss some formal aspects of such a "phase-space density functional theory" and 
explain why it could ameliorate the problems in both cases. For each problem, a two-electron 
model system is exactly numerically solved and analysed in phase-space via the Wigner function 
distribution. 

PACS numbers: 



I. INTRODUCTION 

Time-dependent density functional theory (TDDFT) 
is a remarkably successful theory of many-body systems 
in time-dependent external potentials [1-3]. The analog 
of static density functional theory [4] , TDDFT is based 
on the Runge-Gross proof of a one-to-one mapping be- 
tween the external potential and the time-dependent den- 
sity of electrons evolving under it, for a specified initial 
wavefunction [1] . Because knowledge of the external po- 
tential specifies the many-body Hamiltonian, all proper- 
ties of the interacting electronic system can be extracted 
from just its one-body density and the initial-state. In 
practise, TDDFT utilizes the Kohn-Sham (KS) scheme, 
where the time-dependent density of an interacting sys- 
tem evolving in an external potential Uoxt(r,i) is calcu- 
lated from a fictitious noninteracting system of fermions 
moving in an effective potential, the KS potential Ws(r, t)^ 
that is defined to reproduce the density of the interacting 
system. The KS potential is written as the sum: 



finding the N single-particle orbitals <j)i{r,t), that ini- 
tially made up the initial KS determinant $0. (Atomic 
units, 6^ = 11 — rrie — 1 are used throughout this paper). 
The density of the interacting system is obtained from 
n{r,t) = Runge-Gross theorem 

then assures us that all properties of the true system can 
be extracted from the orbitals. In the linear response 
regime, in which lie the majority of TDDFT applica- 
tions so far, excitation energies and oscillator strengths 
are obtained from a perturbative formulation of these 
equations [5, 6], leading to the matrix formulation that 
operates in most of the quantum chemistry codes today. 

Clearly the accuracy of TDDFT then depends on the 
approximation used for the xc kernel. Despite being 
known generally to have "memory-dependence" - that 
is, Wxc(r,t) depends on the history n{r,t'),t' < t and 
on the initial-states, and $0 [7-9] - almost all ap- 
proximations used today neglect this, and bootstrap a 
ground-state functional: that is, the instantaneous den- 
sity is input into an approximate ground-state functional. 



Vs[n;^o](r,t) = Wcxt(r, t)+WH[n](r, i)+i'xc["; ^o, *o](r, <) w^cl"! *o, *o](r, = v^l[n(rt)](rt) . The superscript A 



where n(r, t) denotes the time-dependent density, and 
\['o('&o) is the initial interacting(nomnteracting) wave- 
function. The second term, Vii[n]{r,t) = / | r-r'| '^^'^^ 
is the classical Hartree potential, while the third is the 
exchange-correlation (xc) potential v^c ■ This is unknown 
as a functional of the time-dependent density and initial- 
states, and must be approximated in practise. Armed 
with an approximation for z)xc['^; *I*o](r, t), one then 
propagates the time-dependent KS equation: 

idtMT^,t) = {-V^/2 + Vs{r,t))(j>,{r,t) 



denotes this "adiabatic approximation" . In the linear re- 
sponse regime [5, 6], adiabatic TDDFT has achieved an 
unprecedented balance between accuracy and efficiency 
for the calculations of excitations and response proper- 
ties. For the first time, one can compute the quantum 
mechanical spectrum of systems as large as biomolecules 
(see eg. Refs. [10-12]), run coupled electron-ion dynamics 
on interesting chemical reactions [13], and study molec- 
ular transport through nanostructures [14]. Excitations 
for which the adiabatic approximation fails have been 
recognized, and, in various stages of being corrected [15- 
(1) 18]. Sometimes the spatial long-ranged- ness of the true 
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functional is important, and approximations more sophis- 
ticated than the usual semi-local ones are required [19- 
22]. 

But TDDFT applies also for strong field dynamics, 

and this is the regime it is particularly promising for, 
given that correlated wavefunction calculations for more 
than two or three electrons in strong-fields become pro- 
hibitively expensive [23]. Compared to calculations of 
spectra, in addition to an approximation for the xc poten- 
tial, now a new ingredient needs to be considered in the 
calculation: observable functionals, i.e. the observables 
of interest must be expressed as functionals of the KS 
wavefunction. The RG theorem guarantees that all prop- 
erties of the true system can be obtained from the KS 
orbitals, as they themselves are implicit density-initial- 
state-functionals, but how? If the property is related 
directly to the density, then no additional observable 
functional is needed. For example, in high-harmonic- 
generation it is the dipole moment, d = / n{r,t)rcfir 
that is of interest, and indeed TDDFT calculations have 
been successful for high-harmonic generation spectra of 
a range of interesting systems [24-26]. But if we are 
interested in, for example, measuring double-ionization 
probabilities, we require the interacting pair-density as 
a functional of the density (or KS orbitals), and this is 
highly non-trivial [27]. Non-sequential double-ionization 
is a fascinating problem that dogged TDDFT for several 
years [27-29]: TDDFT's promise in capturing electron 
correlation made it attractive for this problem, but it 
was soon realized that not only does one need to go be- 
yond the usual semilocal GGA's for the xc functional in 
order to obtain the knee structure [29], but using an un- 
correlated KS expression for the pair-density yielded a 
knee that was too high. (In Ref. [27], an adiabatic cor- 
related expression, based on the ground-state Perdew- 
Wang pair density [30] was shown to lower the knee ap- 
propriately). Approximating a general observable simply 
by the appropriate operator evaluated on the KS wave- 
function can lead to gross inaccuracies, due to the lack 
of correlation in the KS wavefunction itself. It is impor- 
tant to remember that the KS wavefunction is designed 
to reproduce the exact interacting one-body density, but 
is not supposed to be an approximation to the true wave- 
function. Indeed another aspect of the same double- 
ionization problem recently demonstrated this: Wilken 
and Bauer's calculations of ion-momentum-recoil distri- 
butions in double-ionization [31] showed that the KS 
momentum-distributions were drastically wrong, display- 
ing a single maximum instead of the characteristic two- 
hump structure, and with a significantly overestimated 
magnitude. Typically KS momentum distributions are 
not those of the true system, even if the exact KS or- 
bitals were used, that is if the exact Vxc was known for 
the problem, and used to generate the KS orbitals. An- 
other example of this will be given in Sec. Ill A. 

Returning to the question of the approximate 
Vxc [j^; *o , ^o] (r, i) : the accuracy of the adiabatic approx- 
imation for real-time dynamics gets mixed reviews. In 



strong-field double-ionization, the xc potential appears 

not to be significantly non-local in time in a wide range 
of cases, although it depends on how the field is ramped 
on [32]. Dynamics of two electrons in parabolic wells, 
on the other hand, yield an exact correlation potential 
that appear to depend strongly on the history [33 35]. 
We shall revisit Ref. [35] in Sec. IIIB. Although signifi- 
cant advances have been made in understanding [7, 8] and 
modelling [34, 36-40] memory-dependence, it remains to- 
day a difficult problem. 

One situation in which the usual approximations (spa- 
tially semi-local and adiabatic density dependence) do 
poorly, is when the true wavefunction fundamentally can- 
not be described by a single Slater-determinant. This 
is also true for the ground state where molecular dis- 
sociation curves are a notorious problem [41 43]: when 
an electron-pair bond dissociates, the true wavefunction 
develops a Heitler-London character, which in a mini- 
mal description requires two Slater-determinants. The 
KS description however operates via a single Slater- 
determinant. The exact xc potential allows for the fun- 
damental difference in the KS and true wavefunctions 
by developing rather stark and unusual peak and step 
features, that are difficult to capture in approximations. 
In time-dependent problems, coefficients of the interact- 
ing wavefunction expanded in a basis of single-Slater- 
determinants may change in time dramatically: whereas 
at one time a single-Slater determinant may dominate, 
at a later time, two single-Slater determinants are essen- 
tial. Such situations may happen in electronic quantum 
control problems (now becoming experimentally accessi- 
ble with the advent of attosecond lasers) [44]. Consider 
for example evolving the ground-state of the He atom 
(Is'^) to its first accessible excited state (ls2p). The KS 
ground state is a single-slater determinant composed of 
a doubly-occupied spatial orbital. This evolves under 
the KS Hamiltonian, which is a one-body evolution op- 
erator. The target single excitation is however a double 
Slater determinant, and therefore cannot be reached by 
any one-body operator. The exact KS system attains the 
density of the true target state using a single orbital, not 
two. This yields a very unnatural description of the true 
system; consequently, the KS potential develops rather 
unnatural features, difficult to model [8, 9]. The inabil- 
ity of KS methods to change occupation numbers leads 
to challenges in developing approximations. 

The two ingredients needed in a TDDFT calculation 
- approximations for the xc potential and the relevant 
observable functional - determine the accuracy of the re- 
sults. As discussed above, sometimes the simplest ap- 
proximations work well, but in other cases, it is more 
challenging to derive suitable approximations. In this 
paper, we explore a generalization of TDDFT based on 
the phase-space density. That is, instead of taking the 
coordinate-space density, n(r, t) as our basic variable, we 
consider taking the density in phase-space, w{r,p,t): a 
quasi-probability distribution which is somewhat related 
to the probability of finding an electron at position r 
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with momentum p (see Sec. II). The idea here is two- 
fold: First, as the basic variable has more information, 
functionals of them may be simpler to approximate, or 
equally, simple functionals of them may be more accurate 
than simple functionals of the coordinate-density alone, 
while retaining the favorable system-size scaling prop- 
erty of density functional theories. Second, more obscrv- 
ables are directly obtained without needing additional 
observable- functionals, in particular, those pertaining to 
one-body operators, such as momentum, kinetic energy. 
As the coordinate and momentum operators do not com- 
mute, there is no unique phase-space density function, 
and we begin in Sec. II by discussing the Wigner func- 
tion, which turns out to be a particularly useful choice 
for our purposes. We discuss the 1-1 mapping between 
Wigner functions and potentials, for a given initial-state. 
We discuss also the formal equivalence of the Wigner 
function with time-dependent density-matrix functional 
theory, which has seen a recent resurgence of interest [45- 
47]. Next, we return to TDDFT and explore the true 
and KS phase-spaces, and the implications of a phase- 
space based theory, for two examples pertaining to the 
challenges in TDDFT described above: momentum dis- 
tributions (Sec. Ill A) and history-dependence in the xc 
functional (Sec IIIB). 

II. THE 1-BODY WIGNER PHASE-SPACE 
DENSITY 

Quantum phase-space distributions have played an im- 
portant role in the development of cjuantum mechanics. 
Although a precise joint probability distribution in po- 
sition and momentum is impossible due to the Heisen- 
berg uncertainty principle, various quasi-probability dis- 
tributions have been defined, which contain position and 
momentum distribution information in a way consistent 
with the uncertainty principle. Such distributions are 
particularly useful in a semiclassical context, relating 
the quantum states to the underlying classical trajec- 

I 



Wk{ri..rk,pi..pk,t) = 



where the k-th order density ma- 
trix is pfe(ri..rj^,ri..r2;i) = 

(^^^ J^*{r[..r'j^,t)<^{ri..rN,t)drk+i..drN. We 

consider here a density-functional theory based on 
the one-body Wigner function w{r,p,t) = Wi{r,p,t) 



tories, and have been exploited extensively, for exam- 
ple in quantum chaos and quantum optics. Phase-space 
approaches have however been largely, although not en- 
tirely [48-53] neglected in density functional theories, 
where position plays a preferred role over momentum. 
Most recently [51-53], Gill and co-workers, have devel- 
oped and tested models for the ground-state correlation 
energy based on various "phase-space intracules" which 
are essentially different contractions of the two-body re- 
duced Wigner function. For example, there is the Omega 
intracule, which is a function of three variables; the dis- 
tance between two electrons, the magnitude of the differ- 
ence in their momentum vectors, and the angle between 
their position and momentum [53]. This is motivated by 
the very physical idea that the correlation between two 
electrons depends on their relative momentum as well as 
on their relation distance [54]. 

Wigner functions were first introduced in 1932 [55- 
57], with the intent of application to many-body sys- 
tems. The Wigner function quasi-probability distribu- 
tions are real but not positive semi-definite; indeed the 
negative areas have been interpreted as an indication of 
non-classical behaviour, from interference to tunneling 
processes, and research in how to classify these areas con- 
tinues today (see for example Refs. [58]). Although ini- 
tially introduced to treat many-particle dynamics, the 
majority of applications of the Wigner function actu- 
ally involve one-particle systems, especially for making 
quantum-classical correspondence. It has been argued 
that the Wigner function is a particularly suitable ap- 
proach to studying transport [59-61] due to its setting 
within a semiclassical picture while yet being rigorously 
quantum mechanical. Rcf. [61] develops a formalism 
based on reduced- Wigner-functions for indistinguishable 
fermionic systems accounting for Pauli-exchange (but not 
correlation). 

For a system of N identical particles, the reduced k- 
body Wigner phase space density is defined as the Fourier 
transform of the /cth-order density matrix p^: 



I 

(henceforth referred to as simply Wigner function): 



w{r,p,t) = (^^^ J d3ypi(r-y/2,r+y/2,t)e^P-f (3) 



{h) J '^'yi-^'^/fce'^' - yi/2..rk - yk/2, ri + yi/2..rk + yk/2; t) 



(2) 



The one-body densities in coordinate-space and 
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momentum-space are directly obtained by integration: 

n{r,t) = j w{r,p,t)d^p (4) 

h{p,t) = j w{r,p,t)d^r (5) 

Several observations can immediately be made about 
such a phase-space density-functional theory (PSDFT) 
whose basic variable is the Wigner phase-space density 
(taking the role of the coordinate-space density in usual 
density functional theories). First, the theory is formally 
equivalent to the recently developed time-dependent one- 
body density-matrix functional theory (IDMFT) [45-47]: 
w and pi contain exactly the same information, in dif- 
ferent forms related via a Fourier transform. Second, 
as in the IDMFT, there is no Kohn-Sham counterpart: 
the Wigner-function of an interacting system cannot be 
reproduced by any non-interacting system in a local po- 
tential. In a non-interacting system, the one-body re- 
duced density matrix is idempotent, while in an inter- 
acting system, it cannot be. It is interesting that gen- 
erally while a non-interacting system may be found that 
reproduces the coordinate-space density of an interact- 
ing system, and a different non-interacting system may 
be found reproducing its momentum-space density, one 
cannot find a non-interacting system that reproduces its 
phase-space density. Third, the Wigner function directly 
gives the expectation value of any one-body operator: 
no additional observable functional would be needed, for 
example, for kinetic energy or momentum distributions. 
Fourth, the one-to-one mapping between the potential 
and the phase-space density, for a given initial-state, fol- 
lows directly from the Runge-Gross theorem: We have 

w{r,p,t) ^ n{r,t) (6) 
(n(r,t),*o)^*(ri..rjv,i) (7) 

where both arrows indicate a unique mapping, the first 
following from p-integration (Eq. 5), and the second from 
RG, while also 

^{Ti..rN,t) ^w{r,p,t) (8) 

uniquely for a specified initial-state. Hence, for local ex- 
ternal potentials, there is a 1-1 mapping between the one- 
body Wigner function and the many-body wavefunction 
for a given initial state, i.e. that all observables may be 



extracted from the one-body Wigner function and initial- 
state alone. 

We note that this 1-1 mapping holds for local external 
potentials (i.e. purely multiplicative in coordinate-space. 
Gilbert [62] has shown that for the ground state, the ex- 
tension of the Hohenberg-Kohn therem to one-body re- 
duced density matrices applies also to the wider class of 
non-local potentials, v{r,r'). Formally, one may be led 
to search for a one-to-one pi — v{r,r') mapping in the 
time-dependent case also, given that these are conjugate 
variables that couple together in an energy functional, 
(c.f. the n — V mapping of Runge-Gross or Hohenberg- 
Kohn [4] , and the A — j mapping of Ghosh and Dhara in 
TD-current-density functional theory [20, 63]. Although 
a proof for time-dependent non-local potentials is yet to 
be found, the mapping may be generalized to vector po- 
tentials A(r, t), extending the realm of application to in- 
clude magnetic fields: 

w^j{r,t) (9) 
(j(r,t),*o)^*(ri..r;v,i) (10) 

where the first arrow follows from 

j(r,i) = jd^pw{r,p,t) (11) 

and the second from the Ghosh-Dhara proof [63] of the 1- 
1 mapping between currents and vector-potentials. That 
is, any observable of any time-dependent interacting elec- 
tronic system evolving under an external vector poten- 
tial, is a functional of its initial wavefunction and time- 
dependent one-body phase-space density. 

The Runge-Gross proof proceeded in two steps, where 
the first proved a one-to-one mapping between currents 
and (scalar) potentials for a given initial state. The sec- 
ond step, to go from currents to densities, holds imder the 
condition that a boundary term, that involves the density 
and the gradient of the potential, vanishes. The impli- 
cations of this boundary condition have been discussed 
in Ref. [64]; it means that TDDFT cannot be applied to 
periodic systems in uniform electric fields [65]. We point 
out that here there is no requirement on the boundary 
conditions for the w — v (or w — A mappings to hold. 

Time-dependent phase-space functional theory in- 
volves solving the equation of motion for the Wigner 
function: 



w{r,p,t) = 



(p-p')-y 



Vext(r+ |) -t^ext(r- |) 



w{Y,p,t) 



-p • V - i j d^p' j d^ye-' 
- % [ d^yd^pid^p2d^r2e'^P-P'yy (- ^— ^ - , ^"j «;2(r,r2,pi,p2,i) 



(12) 



where the last term must be approximated as a functional of the one-body w and the initial-state ^n- This equation 
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follows directly from that for the density-matrix, which may be more familiar: 
I 

ipi(r', r, t) = (-VV2 + i;ext(r, t) + V'V2 - i^ext(r', t)) pi(r', r,t) + f d'r^ ( . ^ , - yj^ r ) /c>2(r', ra, r, V2,t) 

(13) 

I 



where p2 is the second-order density-matrix, to be ap- 
proximated as a functional of p\. This equation is first in 
the BBGKY (Bogoliubov-Born-Green-Kirkwood-Yvon) 
heirarchy of reduced density matrix evolutions, which 
gives the equation of motion for the fcth reduced DM in 
terms of the [k + l)st and lower DM's. Eq. II is the crux 
of IDMFT, which has only very recently begun to be ex- 
plored. All the development and results so far have been 
in the linear response regime [45-47], using functionals 
for p2 which have been adiabatically bootstrapped from 
ground-state ones (eg. Refs.[66, 67]). 

The accuracy of the time-dependent PSDFT(IDMFT) 
depends on approximate functionals for the term involv- 
ing the second-order Wigner function W2(p2) in terms of 
w{pi). As has been noted above, and earlier in the lit- 
erature, the kinetic term is treated exactly, in contrast 
to DFT, where the correlation potential/energy has a ki- 
netic component. The simplest approximation would be 
the uncorrclatcd one, where p2 is written as a product 
of antisymmctrizcd one-body pi's. This truncates the 
BBGKY heirarchy at the first level, and is equivalent to 
time-dependent Hartree-Fock. Its inability to change oc- 
cupation numbers makes it unable to treat many of the 
challenges that TDDFT faces, that the phase-space DFT 
would hope to treat. For example, the electronic quan- 
tum control problem mentioned in the introduction (see 
also [68]). Section III A discusses another example where 
a change in occupation number is crucial for an accurate 
description. Natural orbitals ipi{r,t) with occupations 
fi{t) are defined by diagonaHzing pi(r, r',t): 

p^{r,r',t)=^Mt)i,i{r,t)ri{r',t) (14) 

i 

We may define "natural Wigner orbitals" via the analo- 
gous expansion: 

w'(r, p,t) =^ fi{t)wi{v, p, t) (15) 

i 

where Wi{Y,Tp,t) = J d^y(j)*{r - y/2)(l)i{r + y/2,t). From 
Eq. II the time evolution of the occupation numbers can 
be derived: 

ifi{t)= j j yp3^P2(r,r2,r',r2,i) 

(i)*{T,t)<i)i{Y' ,t)dTdT' dY-2 - C.C (16) 

The right-hand-side of this equation is only nonzero 
if the contraction of p2 on to the natural orbitals. 



</)*(r, i)(02(r, r2, r', r2, <)</'i(r', f), is imaginary. There- 
fore occupation numbers cannot change with using an 
uncor related product expression for p2- It was also 
shown that bootstrapping any adiabatic functional from 
ground-state IDMFT cannot change occupation num- 
bers [47, 69]. A systematic approach might consider 
then truncating the BBGKY heirarchy at the next or- 
der; that is, solving the equation of motion for p2 con- 
sistently with Eq. II, but using the uncorrelated (anti- 
symmetrized product form) for ps. However, it has been 
shown that this violates fundamental trace relations be- 
tween density matrices of different order [70, 71]; ways to 
get around this become rapidly complicated. Moreover, 
one must deal with a four-point function (second-order 
Wigner function) rather than the two-point one (first- 
order). 

The problem of an appropriate W2-functional in terms 

of w, which is able to change occupation numbers while 
retaining fundamental physical properties, is an impor- 
tant area for future research. 



III. REVISITING TWO CHALLENGES WITHIN 
PHASE-SPACE DFT 

We now describe in a little more detail two problems 

that are particularly challenging for TDDFT approxima- 
tions, and discuss how a PSDFT can ameliorate them. 



A. Kohn-Sham momentum densities 

In general, the momentum-distribution of the KS sys- 
tem is not the same as that of the true system: although 
the sum of the squares of the KS orbitals in coordinate- 
space yield the exact density of the true system, the 
sum-square of their Fourier transforms to momentum- 
space does not yield the exact momentum-density. In 
cases where the correlation component to the kinetic en- 
ergy, Tc is relatively small, calculating the momentum 
distribution simply from the Fourier transforms of the 
relevant KS orbitals can sometimes be quite accurate, 
as recently illustrated by the Compton profiles obtained 
in Ref. [72], provided accurate enough ground-state xc 
potentials are used. Earlier, for the the ground state 
case, Lam and Platzman derived a local-density approx- 
imation for the correlation correction to KS momentum- 
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densities for Compton profiles [73]. The KS orbital mo- 
mentum distribution has also been successful in electron 
momentum spectroscopy [74, 75]: this is essentially a 
triple differential cross-section measurement, where one 
electron is scattered from an atom, molecule, or surface, 
causing the ejection of another electron. The experiment 
obtains a map of the electron momentum distribution be- 
fore ejection, thus imaging the momentum distribution of 
the Dyson orbital. Refs. [74, 75] show that the KS or- 
bital yields a good approximation to the Dyson orbital in 
their momentum distributions, better than, for example, 
Hartree-Fock. 

Away from the ground-state, Wilken and Bauer [31] 
have shown that KS momentum-distributions of TDDFT 
are miserably poor for non-sequential double-ionization 
processes; not only their shape is wrong, lacking the char- 
acteristic double-hump structure, but also their magni- 
tude is significantly underestimated. Ref. [31] develops a 
"product-phase" approximation to extract the true mo- 
mentum distribution from the KS system, i.e. an observ- 
able functional for the momentum-distribution. 

The possibility of working directly in momentum-space 
is hindered by the lack of a one-to-one mapping theo- 
rem for momentum-densities (even for the ground-state 
case) [77, 78]. In fact the earlier interest in developing 
ways to associate phase-space distributions to a given 
ground-state coordinat(^-space density [48-50] was mo- 
tivated by trying to get good momentum-space proper- 
ties derived from the coordinate-space density. One may 
view the density-matrix developments as a rigorous ap- 
proach to this problem; as discussed in Sec. II, the exact 
one-body density matrix, or equivalently, the one-body 
phase-space density, directly yields all one-body informa- 
tion, i.e. the expectation value of any one-body operator, 
such as the momentum-density operator, is exact. 

Even a simple model consisting of two electrons in 
one-dimension illustrates the discrepancy between the KS 
and true momentum-densities in ionization or scattering 
problems. With an exactly solvable model, one can pin- 
point the origin of the differences between the momen- 
tum distributions of the true and KS systems. In this 
particular example, we will see it is due to the single- 
Slater determinant description being a poor description 
of a state in which two electrons are doing completely dif- 
ferent things. One electon is at rest (modelling a "bound" 
electron), while the other is moving away from it (mod- 
elling a "scattered" or "ionized" electron), and we con- 
sider that enough time has elapsed that the true wave- 
function describes independent electrons: 

^{x, x', t) = {tpQ{x)Mx'^ t) + Mx')Mx, t)) 1^ (17) 

The orbital for the bound electron is simply chosen as a 
Gaussian 



i\)q{x) 



(18) 



initially ejected by some laser field (turned off thereafter), 
such that its initial spread in position is 1/A, and average 
momentum is po: 



t/2 j 



(19) 



The coordinate-density, n(x,t) = 2/ dx'\^{x,x',t)\'^, 
and momentum-density, n{p, t) = 

2 / dp'\J dxdx'e*^''+P''''^^{x,x',t)/{2n)\^ of this state 
at a snapshot in time is shown in Fig. 1: Neither are 
remarkable, as both are what might be expected quasi- 
classically from the sum of two probability distributions 
centered around the average position and momentum 
for each electron. Now we consider the Kohn-Sham 
description, which, for two electrons in a singlet state, 
operates via the doubly-occupied orbital, reproducing 
the time-evolving density n{x,t) of the true system: 



(i>{x,t) 



n{x, t) 



exp 



while that for the ionized electron is, for simplicity, cho- 
sen to be a free-particle dispersing Gaussian, which was 



where J{x,t) is the current-density of the true system. 
We compute the KS momentum distribution, ns{p,t) = 
|^(p, i)P from the Fourier transform d of the orbital (f). 
From Figure 1, we sec that while the overall envelope 
of the KS momentum distribution follows that of the 
true system, there are distinct oscillations. In Figure 2 
we plot the corresponding Wigner phase-space distribu- 
tions. w{x,p,t) and Ws{x,p,t); oscillations are evident in 
the latter also, and dip into negative values. This is a 
typical signature of non-classical behavior: the KS or- 
bital describes an electron spatially delocalized in two 
regions, impossible for a classical particle. The coher- 
ence between the two separated parts of the wavepacket 
gives rise to the oscillations in the phase-space density: 
integrating over momentum washes out these oscillations, 
yielding the sum of the two Gaussians in the left-hand- 
side of Fig. 1 (i.e the true interacting density) whereas 
the oscillations persist in the KS momentum distribution 
on the right. On the other hand, the true Wigner func- 
tion looks much like a classical joint phase-space distri- 
bution would be: a peak at the position and momentum 
of each electron, each locally smeared out. Aside from 
the antisymmetrizing, the true wavefunction has classical 
features, with each electron in their independent Gaus- 
sian wavepackcts, quite in contrast to the KS description 
which describes both electrons in the same delocalized 
orbital. Although yielding the exact coordinate density 
of the true system, this exact KS wavefunction is clearly 
far from the true wavefunction. Fundamentally, it asso- 
ciates to the true two-determinant state, a single deter- 
minant; in a sense, this is a time-dependent analog of the 
static correlation problems that haunt ground-state den- 
sity functional theory, for example in dissociation prob- 
lems [16, 41, 76]. We can understand the discrepancy in 
the momentum distributions as arising from this funda- 
mental difference in the nature of the true and KS states 
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in this particular case. (In the next section, we shall see 
an example where the momentum distributions are not 
so different). We stress that this is the KS system with 
the exact xc functional: what is required in TDDFT to 
get the correct exact momentum distribution from this, 
is the (unknown) momentum observable functional. 
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FIG. 1: Density profiles at time f = 8, for parameters A = 0.8 
and po = 2. The left-hand panel shows the coordinate-density 
profile, identical by definition for KS and true. The right-hand 
panel shows the true (solid) and KS (dashed) momentum- 
density profile. 




■0.5 




FIG. 3: The KS Wigner profile at time t = 8, for parameters 
A = 0.8 and po = 2. 



not be there in such a state which is of classical nature. 
We comment here that a similar problem occurs for cer- 
tain quantum control theory type problems, as discussed 
in the introduction. 



B. Dependence on the history of the density 




FIG. 2: The true Wigner profile at time t = 8, for parameters 
A = 0.8 and po = 2. 

What are the impHcations of this for a PSDFT? As 
explained in the earlier section, in a phase-space DFT 
we deal directly with the interacting system; there is no 
non-interacting system that can reproduce an interacting 
phase-space density. Momentum distributions are imme- 
diately obtained without the need for an additional ob- 
servable functional. And, had we the exact functional 
for the W2-term in terms of w, these would be exact. In 
practise, approximations for this term are needed, and 
for this example, it is crucial for the approximation to 
be able to change occupation number. The state starts 
in a ground-state, both electrons bound, which, in a first 
approximation would be in the same doubly-occupied or- 
bital. In the process of ionization, one electron leaves, and 
a two-orbital description is necessary, each with occupa- 
tion number 0.5. Any approximation unable to change 
occupation numbers will once again yield wiggles in the 
momentum and phase-space distribution. These should 



The time-dependent exchange-correlation potential 
functionally depends on the initial state (both the true 
and KS) as well the history of the density [1, 8, 
9]. This dependence on the past is referred to as 
"memory-dependence" . Any adiabatic approximation 
(e.g., ALDA) utilizes only the instantaneous density, 
so ignore this temporal nonlocality of the xc poten- 
tial. Yet, adiabatic approximations work remarkably 
well for the calculation of most, but not all, excitations 
and their oscillator strengths. For dynamics in strong 
fields, there arc both parameter regimes and phenomena 
in which the memory-dependence is minimal, but also 
examples where memory-dependence is vital [8, 9, 33]. 
Memory-dependence is not easy to capture in approx- 
imations although there has been intense progress re- 
cently [34, 36, 38-40]. 

Ref. [35] used a numerically exactly-solvable example 
to demonstrate explicitly a case where the correlation po- 
tential depends on the density ultra-nonlocally in time. 
There, time slices were found in which the density was 
practically identical, semi-locally in time, whereas the 
correlation potential differed significantly. We shall re- 
visit this example here, with a view to asking whether 
functionals in a PSDFT could be less non-local in time, 
and thus easier to approximate. Specifically, we shall ask 
whether the phase-space density can distinguish those 
time slices mentioned above, that the coordinate-space 
density could not. We shall find that indeed it can, and 
moreover, that the size of the differences in the phase- 
space density tracks the difference in the correlation po- 
tential. This suggests that the memory-dependence in a 
PSDFT could be milder than in the usual DFT. 

Time-dependent Hooke's atom is a system of two 
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Coulombically interacting electrons in a time-dependent 
parabolic well: the Hamiltonian is 



H = -Iv? - + \k{t){rl + rl) + —1— (21) 



For our investigations, wc chose k(t) = k — ecos{Lot), as 
in Ref. [35] . The Hamiltonian decouples into relative (w) 
and center of mass {R) co-ordinates , so one may solve 
for the exact time-dependent 2-electron wave function 
as a product of R and u wavefunctions, 5'(ri,r2,t) = 
xiR, t) and thereby obtain time-dependent density 

n{r,t). Knowledge of the time-evolution of the density 
yields the doubly-occupied time-dependent KS orbital, 
and then, via inversion of the time-dependent KS equa- 
tion, the exchange-correlation potential Vxc{r,t). A de- 
tailed derivation can be found in Refs. [35, 79], along with 
computational details. We mention here only that the 
time-propagation method used is the Crank-Nicholson 
scheme, and a grid chosen to sample points near the ori- 
gin more densely than further away, since the density 
decays exponentially. With chosen parameters e = 0.1, 
w = 0.75, and k = 0.25, the density retains a near- 
Gaussian shape at all times, with a time-varying width; 
thus it may be parametrized by its variance rrms{t) = 
a/ (r^) [35]. This is shown in the top panel of the Fig. 4. 
The lower panel shows a density- weighted correlation po- 
tential, Ec = J d^rh{r,t)vc{r,t), which was also studied 
in the earlier paper [35], as a simple way to track the 
time-evolution of the correlation potential. A pair of time 
slices is indicated within which rrms is practically identi- 
cal in each time slice, while the density-weighted correla- 
tion potential Eq differs quite dramatically. Other such 
pairs of time slices may be found. This indicates that 
local, or even semi-local in time density-dependence is 
not enough to specify the correlation potential, i.e. that 
the correlation potential depends on a significant history 
of the density [9, 35]. Any adiabatic approximation that 
uses only instantaneous density information would erro- 
neously predict the same correlation potential for each 
time slice. 

We now consider the phase-space distribution at these 
time slices: if these are also nearly identical at these 
pairs of time slices, then it would suggest the memory- 
dependence in a PSDFT is just as non-local as in DFT 
and just as challenging to model. We shall instead find 
that the phase-space distributions distinguish the system 
at these difi'erent pairs of time slices, with the differences 
greater with greater differences in the density-weighted 
correlation potential Ec ■ This suggests that perhaps the 
memory-dependence would be less severe in PSDFT, and 
therefore perhaps easier to approximate. 

Knowledge of the true and KS wavefunctions, enables 
us to derive expressions for Wigner phase-space distribu- 




FIG. 4: Nonlocality of the correlation potential in time: 

J Vc{r,t)h{r,t)d^r (lower panel) and rrms(t) (top panel). In 
the two time slices indicated, the density is almost identical 
locally in time, whereas the donsity-woightod correlation po- 
tential in the lower panel is quite distinct. 



tion for the interacting and KS cases: 



w{r,p,t) 



/ x^dx I 



-ipx cos a cos 9 



Jo{px sin ay/l- cos^ 6')x*(^v^' *)^(v^' t)d{cos 9) 



(22) 



-ipx cos a cos 



POO /* I 

t«s(r, p, = / x^dx / e" 
Jo J -I 

X Jo{2px smay/ 1 — cos^ 0)0* (iy^, t)4>{y/rZ, t)d{cos fl^3) 

where = r'^ + x'^ + 2rx cos 6, r_ = + — 2rx cos 9 
and a is the angle between r and p and Jq is the Bessel 
function of zeroth order. In Eq. 23, <j) is the KS orbital. 
Both distributions w and Ws depend only on the mag- 
nitudes of the coordinate and momentum, r and p, and 
the phase-space angle a; for simplicity, we look at the 
angle-averaged quantities: 



.j{r,p,t) = / w{r, p, a, t) sin ada (24) 
Jo 



Further, to highlight the component of the Wigner func- 
tion that is due to correlation, we take the difference be- 
tween the angle-averaged true and KS Wigner functions, 
defining 

Wc {r, P, t) = Wav {r, p, t) - Ws,av {r, P, t) (25) 

To simplify the analysis even further, we consider only 
the momentum-distribution of this correlation compo- 
nent, defined as: 

.„(p,<)^*(p,,)-*.(p,,)^/».(.,.,,y* (26, 
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FIG. 5: Top left panel: the true (solid) and KS (dashed) 

momentum densities, at t = 29,8. Top right and lower two 
panels: correlation components of the momentum distribu- 
tions nc{p) at three different time-pairs indicated. In each 
time-pair, the coordinate density, n{r,t), is identical, while 
the momentum-density is not, with the difference in its cor- 
relation component nc(p) growing as the difference in the 
density-weighted correlation potential Ec at those times. 



We now turn to the plots in Figure 5. The top left panel 
shows the true and KS momentum probabiUty distribu- 
tions at t = 29.8, which has one of the larger values 
for the density-weighted correlation potential in this run 
(Fig. 4). Wc include in all these plots the Jacobian fac- 
tor, so that we are actually plotting p'^n{p). Although 
the momentum distributions are not identical, the differ- 
ences are relatively small, in contrast to the ionization 
example in the previous section. We now study their 
difference ndp)- The other three panels show ndp) at 
pairs of times at which the coordinate-density are prac- 
tically identical (i.e. essentially indistinguishable to the 
eye on the scale of plots similar to the top left panel). 
The first time-pair, (4.8,29.8), is chosen near the mini- 
mum of the r„ns at each of the two time slices indicated 
in Figure 4; here the difference in Ec is relatively small 
(0.0005 au), while it is bigger (0.011 au) at the second 
time pair (5.1,29.3), and bigger still (0.0264) at the third 
pair (5.5,29.3). In contrast to the coordinate-density, the 
momentum-densities are not identical at these pairs of 
times. Moreover, their difference increases with the dif- 
ference in the value of the density-weighted correlation 
potential at those times. The momentum distribution 
appears to somewhat track the correlation potential. 

What are the implications of this for memory in a PS- 
DFT? The results suggest that having momentum infor- 
mation in the basic variable may reduce their memory- 
dependence. Of course this is not conclusive, but the fact 
that momentum distributions distinguish the system at 
pairs of times when the coordinate-density does not, yet 
the correlation potential is different, does suggest that 
functional of the phase-space density may be less non- 
local in time. 



IV. SUMMARY AND OUTLOOK 

A PSDFT approach could be more successful than the 
usual TDDFT for some applications for which the usual 

functional approximations in TDDFT perform poorly. 
In many strong-field applications, memory-dependence 
is required in the TDDFT xc potential, challenging to 
model, and additional observable functionals arc needed, 
depending on the measured quantity of interest. The ex- 
amples worked through in section III demonstrate each of 
these aspects clearly. Being exactly-numerically solvable 
two-electron problems, one can perform a detailed anal- 
ysis that helps to understand these challenges in their 
broader context. 

The first example, a model of ionization, illustrated, 
via a phase-space exploration, that even though an ex- 
act KS treatment gets the exact coordinate-space density 
correct, it may do so with a wavefunction that is fun- 
damentally different than the true wavefunction; there- 
fore approximating observables not directly related to 
the coordinate density by their value on the KS or- 
bitals does not work well. Wc considered the momen- 
tum distribution (eg modelling that as measured by an 
ion- momentum recoil experiment), obtained from the KS 
orbitals and found spurious oscillations as a consequence 
of the single-determinantal nature of the KS wavefunc- 
tion: a single delocalized orbital is describing the true 
two-orbital state. The exact coordinate-space density is 
recovered (by definition), but the momentum distribu- 
tion displays oscillations, present also in the KS phase- 
space distribution. The latter goes negative, a sign of 
non-classical behavior: the KS system describes the ex- 
act density putting both electrons in a doubly-occupied 
orbital, which must therefore be delocalized in two re- 
gions. The true system is however a rather classical one, 
aside from the antisymmetry, and this is reflected in its 
phase-space density profile. 

The second example explored the memory-dependence 
of the correlation potential in a time-dependent Hooke's 
atom. Solving directly for (a density-weighted measure 
of) the corrc;lation potential, time slices were found in 
which the coordinate-density evolved identically, while 
the correlation potential varied dramatically, illustrating 
the fact that the correlation potential is a highly non- 
local in time functional of the coordinate-density. We 
then asked whether, if momentum-distribution informa- 
tion was included in the basic variable, the correlation 
potential would be just as non-local in time as a func- 
tional of the phase-space density, or whether the phase- 
space density could distinguish the state of the system 
at these times. We found the latter was true, suggesting 
that correlation functionals of the phase-space density 
are less non-local in time, and therefore easier to approx- 
imate than correlation functionals of the density alone. 

In PSDFT, the basic variable is the one-body Wigner 
phase-space density, from which momentum-space and 
position-space densities can be obtained via integration 
over the conjugate variable. Formally equivalent to 
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density-matrix functional theory, all one-body observ- 

ablcs arc directly obtained from the basic variable (in- 
cluding the momentum distribution), without the need 
for additional observable functionals. In its equation of 
motion the action of the kinetic energy operator is exact 
as a phase-space functional (equally, density-matrix func- 
tional), while approximations are needed for the second- 
order Wigner function as a functional of the first-order 
one. An important future direction is to develop func- 



tional approximations for this that lead to changing oc- 
cupation numbers, in order to treat many of the prob- 
lems that are challenging in the usual TDDFT, such as 
electronic quantum control problems, and ionization dy- 
namics. In light of the studies in this paper, this would 
be well worth pursuing. 
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